/*
 * DistanceFunction.cpp
 *
 *  Created on: Jan 15, 2011
 *      Author: LyonsDesktop
 */

#include "DistanceFunction.h"
#include "../../Datasets/Dataset.h"

#define USE_SSE

#ifdef USE_SSE
#  include "../../SSE.h"
#endif

//Frechet distance

using namespace CLEVER;
using namespace CLEVER::Datasets;

void CLEVER::RegionDiscovery::DistanceFunctions::MinkowskiDistanceFunction::DistanceSSE(NumericDataset * d1, int i1, int i1b, DoubleArrayDataset * d2, int i2, double & di12, double & di13, bool Relative)
{
	//	std::cout<<"In DistanceSSE(DoubleArray)"<<std::endl;

	double udist1 = 0;
	double udist2 = 0;

	int num_attr = d1->GetNumAttributes();

	double * p1 = d1->GetValuePtr(i1);
	double * p1b = d1->GetValuePtr(i1b);
	double * p2 = d2->GetValuePtr(i2);

	/*
	   #ifdef USE_SSE
	    if (num_attr == 2)
	    {
	        __m128 expw = _mm_set_ps1((float)this->exp);
	        __m128 invexpw = _mm_set_ps1((float)invexp);

	        __m128 a = _mm_set_ps((float)p1b[1], (float)p1b[0], (float)p1[1], (float)p1[0]);
	        __m128 b = _mm_set_ps((float)p2[1], (float)p2[0], (float)p2[1], (float)p2[0]);
	        a = _mm_abs_ps(_mm_sub_ps(a, b));

	        __m128 v = _mm_exp_ps(_mm_mul_ps(expw, _mm_log_ps(a)));
	        v = _mm_hadd_ps(v, v);
	        if (!Relative)
	        {
	            v =  _mm_exp_ps(_mm_mul_ps(invexpw, _mm_log_ps(v)));
	        }

	        __declspec(align(16))
	        float f[4] = { 0 };
	        _mm_store_ps(f, v);

	        di12 = f[0];
	        di13 = f[1];
	        return;
	    }
	   #endif
	 */

	for (int i = 0; i < num_attr; i++) {
		//double x = d1->GetValue(i1, i);
		//double y = d2->GetValue(i2, i);
		double diff = fabs(p1[i] - p2[i]);
		udist1 += pow(diff, exp);
	}
	for (int i = 0; i < num_attr; i++) {
		//double x = d1->GetValue(i1, i);
		//double y = d2->GetValue(i2, i);
		double diff = fabs(p1b[i] - p2[i]);
		udist2 += pow(diff, exp);
	}
	di12 = pow(udist1, invexp);
	di13 = pow(udist2, invexp);
}


void CLEVER::RegionDiscovery::DistanceFunctions::MinkowskiDistanceFunction::DistanceSSE(NumericDataset * d1, int i1, int i1b, NumericDataset * d2, int i2, double & di12, double & di13, bool Relative)
{
	double udist1 = 0;
	double udist2 = 0;

	int num_attr = d1->GetNumAttributes();

	double * p1 = d1->GetValuePtr(i1);
	double * p1b = d1->GetValuePtr(i1b);
	double * p2 = d2->GetValuePtr(i2);

	/*
	   #ifdef USE_SSE
	    if (num_attr == 2)
	    {
	        __m128 expw = _mm_set_ps1((float)this->exp);
	        __m128 invexpw = _mm_set_ps1((float)invexp);

	        __m128 a = _mm_set_ps((float)p1b[1], (float)p1b[0], (float)p1[1], (float)p1[0]);
	        __m128 b = _mm_set_ps((float)p2[1], (float)p2[0], (float)p2[1], (float)p2[0]);
	        a = _mm_abs_ps(_mm_sub_ps(a, b));

	        __m128 v = _mm_exp_ps(_mm_mul_ps(expw, _mm_log_ps(a)));
	        v = _mm_hadd_ps(v, v);
	        if (!Relative)
	        {
	            v =  _mm_exp_ps(_mm_mul_ps(invexpw, _mm_log_ps(v)));
	        }

	        __declspec(align(16))
	        float f[4] = { 0 };
	        _mm_store_ps(f, v);

	        di12 = f[0];
	        di13 = f[1];
	        return;
	    }
	   #endif
	 */

	for (int i = 0; i < num_attr; i++) {
		//double x = d1->GetValue(i1, i);
		//double y = d2->GetValue(i2, i);
		double diff = fabs(p1[i] - p2[i]);
		udist1 += pow(diff, exp);
	}
	for (int i = 0; i < num_attr; i++) {
		//double x = d1->GetValue(i1, i);
		//double y = d2->GetValue(i2, i);
		double diff = fabs(p1b[i] - p2[i]);
		udist2 += pow(diff, exp);
	}
	di12 = pow(udist1, invexp);
	di13 = pow(udist2, invexp);
}


double CLEVER::RegionDiscovery::DistanceFunctions::MinkowskiDistanceFunction::Distance(NumericDataset * d1, int i1, NumericDataset * d2, int i2)
{
	double udist = 0;
	int num_attr = d1->GetNumAttributes();

	double * p1 = d1->GetValuePtr(i1);
	double * p2 = d2->GetValuePtr(i2);

	//#ifdef USE_SSE
	//	if (num_attr == 2)
	//	{
	//		__declspec(align(16))
	//		float f[4];
	//
	//		__m128 a = _mm_set_ps(1, p1[1], p1[0], 1);
	//		__m128 b = _mm_set_ps(0, p2[1], p2[0], 0);
	//		a = _mm_abs_ps(_mm_sub_ps(a, b));
	//		b = _mm_set_ps(exp, exp, exp, exp);
	//		__m128 v = exp_ps(_mm_mul_ps(b, log_ps(a)));
	//		_mm_store_ps(f, v);
	//
	//		return powf(f[1] + f[2], invexp);
	//	}
	//#endif

	for (int i = 0; i < num_attr; i++) {
		//double x = d1->GetValue(i1, i);
		//double y = d2->GetValue(i2, i);
		double diff = fabs(p1[i] - p2[i]);

		udist += pow(diff, exp);
	}
	return pow(udist, invexp);
}


